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Abstract: We investigate the dynamics of vortex formation and 
the structure of vortex lattices in a Bose-Einstein condensate con- 
fined within a rotating double-well (DW) potential. The dynam- 
ical process is a formation process of "ghost" vortices, "hidden" 
vortices and "visible" vortices. The critical rotation frequency for 
the creation of visible vortex is indicated by a sudden jump in the 
angular momentum evolution and an inflexion in the energy evo- 
lution. Different visible vortex structures can be formed by ruling 
the rotation frequency. In particular, structural change of visible 
vortex patterns can be achieved by regulating the anisotropy pa- 
rameter of the DW potential. This feature allows to flexibly con- 
trol the distribution of angular momentum in macroscopic quan- 
tum systems and study the interplay among rotation, interparticle 
interaction and external potential in superfluids. 
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Linear visible vortex lattice pair in a Bose-Einstein condensate 
trapped in a rotating double-well potential 
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1. Introduction 

Topological defects have been explored extensively in 
many systems, such as superfluids 4 He HI and 3 He 13, su- 
perconductors [3], and cosmology [4|. Recently, increas- 
ing interest has been focused on the topological defects 
in ultra-cold atomic gases [.5,6,7]. As an important exam- 
ple of topological defect, quantized vortex plays a funda- 
mental role in the studies of atomic gases, ranging from 
manifesting superfluidity [8], collective excitations, long- 
range phase coherence fljJ). and universal phase transition 
dynamics of Bose-Einstein condensates (BECs) and de- 
generate Fermi gases ifTUI . In particular, for the vortices in 
dilute BECs, the theoretical description can be carried out 
in a much more efficient way than that in liquid He due to 
the low density and weak interatomic interaction. 

To generate vortices in trapped BECs, different ap- 
proaches, such as stirring a BEC with a laser HI 111121 . 
phase imprinting, rotating traps ifPTl . and vector gauge po- 
tential, have been used. Here we concentrate on the case 



of rotating traps or equivalently rotating BECs. It is well 
known that for a BEC in a rotating harmonic trap there is 
a critical rotation frequency beyond which vortex starts to 
exist 1:6] . With the increasing of rotation frequency, more 
and more vortices appear. Structure of vortices depends 
on external trap potential. The vortices form a triangular 
lattice in a symmetry harmonic trap, i.e., the celebrated 
Abrikosov lattice ifTZl . While for rotating anisotropic har- 
monic traps [14.1 15.1 16,| or atomic waveguides ifTTl . linear 
vortex lattices 010, zigzag arrangements lfl9l and stripe- 
shaped configurations 1 20 21 22 1 might be generated. Fur- 
thermore, the vortex patterns undergo structural changes 
depending on the anisotropic parameter. Square vortex lat- 
tices and structural phase transitions (from triangular lat- 
tices to square lattices) has also been investigated in the 
system consisting of harmonic trap and optical lattice 11231 . 
Meanwhile, the combined harmonic and quartic traps sup- 
port giant vortices 11241 . The toroidal traps can hold persis- 
tent flows with high winding numbers l25l . 
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To the best of our knowledge, most of the existing 
theories in the literature mentioned above directly use an 
imaginary time propagation (ITP) method or a lowest Lan- 
dau level (LLL) approximation or a variational approach 
to study the properties of a rotating BEC. These methods 
are incapable of revealing the dynamics of vortex forma- 
tion although they can be applied to find the steady states 
of a rotating system. In addition, these approaches can 
not deal with the BECs with dissipation that is universal 
and inevitable in real cold-atom experiments. Recently, a 
phenomenological dissipation model has been presented to 
study the dynamics of vortex formation of a BEC in a ro- 
tating harmonic trap, where elusive "ghost" vortices were 
revealed l26l . In a most recent investigation ll27l . the au- 
thors proposed an improved phenomenological dissipation 
model to test the well-known Feynman rule Il6l l28l of vor- 
tices in a rotating double-well (DW) BEC. It is found that 
there are three fundamental types of vortices in a rotating 
DW BEC: "visible" vortex, "ghost" vortex and "hidden" 
vortex. Only after including the hidden vortices can the 
Feynman rule be satisfied. 

In this Letter, we study the dynamics of vortex for- 
mation in a BEC trapped in a rotating DW potential by 
using an improved phenomenological dissipation model 
E71 . The DW potential is particularly important due to its 
simplicity, universality and yet richness 129113011311 . Mean- 
while, it is a natural and flexible anisotropic trap. We dis- 
cuss the structural change of visible vortex patterns in a 
rotating DW BEC by regulating the rotation frequency and 
the DW configuration. We find that the dynamical process 
of vortex formation in a rotating DW potential is remark- 
ably different from that in a rotating harmonic trap 11261 . 
The former case is a formation process of ghost vortices, 
hidden vortices and visible vortices. The critical rotation 
frequency for the generation of visible vortices is signaled 
by a sudden jump in the evolution of the average angular 
momentum per atom and an inflexion in the evolution of 
the average energy per atom. It is shown that the rotating 
drive excites a complex turbulent oscillation mainly con- 
sisting of the surface mode with I = 4 rather than the usual 
quadrupole oscillation consisting of quadrupole mode with 
I = 2 in the case of rotating harmonic trap ll26l . Due to 
the parity effect, the vortex patterns of the rotating DW 
BEC in an equilibrium state display well centrosymmetry. 
By controlling the rotation frequency and the anisotropic 
parameter of the DW potential, the visible vortex pattern 
undergoes structural changes, including the formation of 
linear and zigzag vortex lattice pairs. In particular, the lin- 
ear visible vortex lattice pairs have potential applications 
in quantum information processing due to the advantages 
of energetic stability, long-time coherence, controllable in- 
teraction and the least number of nearest neighbors. 



2. Phenomenological dissipation model 

By assuming strong confinement in the z direction, we 
consider a two-dimensional (2D) BEC trapped in a DW 



potential 



U(x,y) = '^{u 2 x x 2 +u*y*) + V e 
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-x 2 /2S 2 
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where m, is the atomic mass, ui x and oj v are the angular fre- 
quencies of the harmonic trap in x and y directions, respec- 
tively. We introduce an anisotropy parameter A = Lo y /tu x 
to describe the ratio between ui y and ui x . Uo and 5 denote 
the height and width of the potential barrier, respectively. 
The DW potential is allowed to rotate around the z axis 
with angular frequency Q. In the rotating frame, the en- 
ergy functional of the system becomes 



E W*} 
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where &(x,y,t) is the wave function of the BEC, and 
L z = ih(y-^ — x-^j) denotes the z component of the 

angular- momentum operator. g2D = 2\f2~nh 2 a s / {a z rn) is 
the 2D interatomic interaction strength with a s being the 
s-wave scattering length and a z = y/h/ (mu> z ) the axial 
harmonic length. If the dissipation of the BEC due to ther- 
mal clouds and rotation are ignored, one can obtain the 
time-dependent Gross-Pitaevskii (GP) equation in the ro- 
tating frame by using a variational method as follow 



lh ~dt = 



~^{Vl + V 2 y )+U(x,y)-QL z 



(3) 



The wave function is normalized as N = J J \>P\ 2 dxdy 
with N being atom number. The stationary solutions 
W(x,y,t) = <P(x,y)e~' l ^ t t h satisfy the time-independent 
GP equation 



fi<P 



2m 



(V 2 x + V 2 v )+U(x,y)-f2L 2 



-92dW\ $, 



(4) 



where fj, is chemical potential. 

Equation (0 or the imaginary time propagation of the 
time-dependent version (0 are widely used to study the 
properties of the stationary states of a rotating BEC. How- 
ever, Eq. (0 and Eq. d4j can not reveal the dynamical pro- 
cess of the vortex formation. In Ref. 11321 . the authors nu- 
merically solved the time-dependent GP equation in a ro- 
tating frame, but the motion of created vortices remains 
turbulent and no vortex lattice forms. Moreover, Eq. (0 
and Eq. (01 can not describe a BEC with dissipation. 

Here we use a phenomenological dissipation model 
||27|| to study the dynamics of vortex generation and the 
structure of the equilibrium state of rotating DW BEC. In 
this phenomenological model, the wave function W(x,y,t) 
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obeys the following modified GP equation 
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where 7 describes the degree of dissipation of the BEC. 
This model is a variation of that in Ref. [26 1 and has good 
predictive power. Our computation results for the case of 
a BEC in rotating harmonic trap are well in agreement 
with the experimental observations in Ref. iTPTl and the 
numerical results in Ref. [26 1 . Recently several other phe- 
nomenological models have also been proposed to deal 
with the BECs with dissipation H33II34U35U36H . To numeri- 
cally simulate Eq. (f5J), it is convenient to introduce dimen- 
sionless parameters as following do = yjh/2mui x , x = 
x/d ,y =jj/d ,cr = SJd ,t = tu> x ,V = U/fkj x ,V = 
U /hu x ,Q = {2/lu x ,L z = L z /h and ip = Wdo/y/N. 
Thus Eq. (|5]l is reduced to a dimensionless form as 

dtp 
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at 
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V(x,y)-f2L; 



where c 



+c\ip\ V, 
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4) 



denotes the dimensionless 
nonlinear strength and the tilde is omitted for simplicity. 
The rescaled 2D DW potential is expressed by 



V(x,y) = ±(x 2 + \ 2 y 2 ) + V e-* 2 / 2 ° 



(7) 



8 ' Rb atoms and the parameters are 
2tt x 40 Hz, uj z = 2tt x 800 Hz, 
Here we take 7 



We numerically solve Eq. © with the split-step 
Fourier method [37.38]. The initial quantum state 
ijj(x,y,t = 0) in the DW potential can be obtained by 
the ITP method I39,40j411 for Q = 0. In our numerical 
computations, we use 
chosen as oj x = u> y = 

Vq = 40 and a = 0.5. Here we take 7 = 0.03, which 
corresponds to a temperature of around 0.1T C [33). In fact, 
the variation of the nonzero 7 in Eq. © only influences the 
relaxation time scale but does not change the dynamics of 
the vortex formation and the ultimate steady structure of 
the system. 



3. Dynamics of vortex formation 

To reveal the dynamics of vortex formation in the rotat- 
ing DW BEC for a fixed !?, we show the time evolution 
of the density distribution \ tp\ 2 and that of the correspond- 
ing phase distribution of tp(x, y, t) in Fig. 1 and Fig. 2, re- 
spectively. In Fig. 2, the value of the phase varies contin- 
uously from to 2-7T, and the end point of the boundary 
between a 2ir phase line and a phase line corresponds 
to a phase defect (i.e., a vortex). Initially, the condensate 
density \ip(x, y,t = 0)| 2 in a stationary DW potential is 
displayed in Fig. 1(a), and there is a zero relative phase 
between the two halves of the BEC as shown in Fig. 2(a). 



With the development of time, the BEC undergoes com- 
plex turbulent oscillation [see Fig. 1(b)] instead of the 
usual elliptic oscillation in the case of rotating harmonic 
trap 11261 . In Ref. l26l . the elliptic oscillation is mainly 
caused by the quadrupole mode with angular momentum 
I = 2. Since the DW system has even parity with respect 
to the spatial coordinate, no surface modes with odd I are 
excited. Here the turbulent oscillation principally consists 
of the surface mode with I = 4 coupling with higher- 
energy modes through the nonlinear interaction, which is 
verified by the critical rotation frequency for visible vor- 
tex generation. The turbulent oscillation makes the bound- 
ary surfaces of the BEC unstable and excites the surface 
waves that propagate along the surfaces. Essentially, the 
surface ripples are resulted from the dynamical instability 
ll42l associated with the imaginary frequency of the excita- 
tion modes. Then the surface waves develop into the vortex 
cores and most of the phase defects appear at the boundary 
surfaces of the cloud as shown in Fig. 2(b). Since the am- 
plitude |-0| on the outskirts of the BEC is almost negligible, 
these phase defects can not be seen in the density distribu- 
tion of Fig. 1(b) and neither carry angular momentum nor 
energy, they are called ghost vortices 112611271 . 

At the same time, there are some phase singulari- 
ties congregating and distributing symmetrically along the 
central barrier [see Figs. 2(b)-2(c)]. These phase singular- 
ities are known as the hidden vortices [27l| because they 
carry evident angular momentum despite being invisible 
in the density distributions of Figs. 1(b)- 1(c). With the fur- 
ther time evolution some ghost vortices penetrate into the 
BEC via the Landau instability [43] associated with the 
negative excitation frequency and becomes the usual visi- 
ble vortices as shown in Figs. 1(d)- 1(e) and Figs. 2(d)-2(e). 
In the presence of dissipation, a steady triangular lattice 
pair of visible vortices [see Fig. 1(f)] forms eventually due 
to the competition between the rotating drive propelling 
vortices toward the rotation axis and the repulsive interac- 
tion tending to push the vortices apart, where the energy of 
the system reaches the minimum in the rotating frame. The 
visible vortex density is determined by the parameters fl, 
A, Vq, <j and c. As displayed in Fig. 2(f), most of the ghost 
vortices are repelled to the outside of the cloud. Since the 
DW potential is even parity concerning to the spatial co- 
ordinate, all the numbers of the visible vortices, hidden 
vortices and ghost vortices should be even, which can be 
seen clearly in Fig. 1(f) and Fig. 2(f). Physically, since vor- 
tex is a topological singularity, the system undergoes the 
change of topological property between the BEC with vor- 
tices and that without vortices when the external rotating 
drive inputs angular momentum into the BEC. Therefore 
the vortices are easier generated in the area of the lower 
density, such as the barrier boundary and the outskirts of 
the cloud [see Figs. 1(b)- 1(c) and Figs. 2(b)-2(c)]. This is 
the reason why the hidden vortices and ghost vortices al- 
ways appear before the visible vortices. 

In Figs. 3(a)-3(b), we display the time evolution 
of the average angular momentum per atom l z = 
J J ip* L z ipdxdy and the average energy per atom e = 
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Figure 1 (online color) Time evolution of the density distribu- 
tion |^| 2 after the double- well potential with A = 1 suddenly 
begins to rotate with Q = 0.9, where c = 600. The time is (a) 
t = 0, (b) t = 10, (c) t = 50, (d) t = 100, (e) t = 160, and (f) 
t = 250. The darker color area indicates the lower density. Here 
x and y are in units of do, and t is in units of 1 /u> x . 
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Figure 2 (online color) Time evolution of the phase distribution 
of ip(x, y,t) corresponding to that of the density distribution of 
Fig. 1. The value of the phase varies continuously from to 2tt. 
The darker color area indicates the lower phase. 



// r h(V= + Vj) + V(x, y) + (c/2) - QL^dxdy, 
respectively. Here A = 1, and Q or c have different 
values. When Q < 0.59, we can see that l z under- 
goes damped oscillation until a convergent value and e 
decreases monotonously to an equilibrium value. When 
Q > 0.59, however, l z first experiences damped oscilla- 
tion, then increases to a new peak value, and finally de- 
scends damply to an equilibrium value [see Fig. 3(a)]. At 
the same time, e first decreases, then passes over an en- 
ergy barrier, and finally drops gradually to a convergent 
value [see Fig. 3(b)]. This distinct difference indicates that 
the critical rotation frequency for the creation of visible 
vortices is fl c = 0.59 which characterizes the topological 
phase transition of the rotating DW system. The critical 
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Figure 3 Time evolution of (a) the average angular momentum 
per atom l z and (b) the average energy per atom £ for A = 1. 
Here t, l z and e are in units of l/to x , ft and huj x , respectively. 



frequency fi c = 0.59 agrees well with that predicted in 
Ref. | [27ll . However, this value is lower than the usual crit- 
ical frequency Qhc = \/2/2 — 0.707 predicted for a ro- 
tating harmonic trap due to the dynamical instability [26. 
l44l . In view of the even parity of the DW potential, the 
most probable reason is that here the rotating drive excites 
mainly the high-order surface mode with I = 4 rather than 
the quadrupole mode with I = 2, where Qj lc = \fljl. The 
small difference between Q c = 0.59 and Qhc = 0.5 is due 
to the nonlinear interatomic interaction and the presence 
of the Gaussian barrier as well as the resulting equivalent 
anisotropy. Moreover, the larger is the Q (or c) with the 
same c (or Q), the larger is the l z in the equilibrium state 
and the lower (or higher) is the corresponding e with a 
slower (or faster) decline, which reveals that the vortex ar- 
ray is a kind of locally stable quantum state. As shown in 
Fig. 1 and Fig. 2, with the development of time the hid- 
den vortices and some ghost vortices penetrate toward the 
areas of high density in the cloud. Due to the DW poten- 
tial consisting of a harmonic trap and a Gaussian barrier, 
the vortex formation process is an oscillating dynamical 
process concerning to the inputting of the external angular 
momentum, which can explain the oscillation behavior of 
the angular momentum in Fig. 3(a). 



4. Structural change of vortex patterns 

The variation of the rotation frequency influences not only 
the dynamics of the vortex formation but also the equilib- 
rium structure especially the visible vortex pattern of the 
rotating DW system. In Fig. 4, we display the steady den- 
sity distributions |')/ , | 2 at t = 250 after the DW potential 
rotates with different angular frequencies, where A = 1 
and c = 600. The critical rotation frequency at which vis- 
ible vortices can enter the BEC with A = 1 is fl c = 0.59 
as has been pointed out above and in Ref. l27l . From Fig. 
4(a), there is a visible vortex pair for the case of fi = 0.65, 
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Figure 4 (online color) Steady density distributions at t — 
250 for different rotation frequencies: (a) Q = 0.65; (b) Q — 
0.7; (c) n = 0.75; (d) Q = 0.8; (e) O = 0.85; (f) Q = 0.95. 
Here A = 1 and c = 600. The darker color area indicates the 
lower density. 



Figure 5 (online color) Steady density distributions \tp\ at t = 
250 for different values of A, where c = 600. (a) A = 2, Q = 
0.9, (b) A = 1.5, n = 0.9, (c) A = 1.1, 1? = 0.9, (d) A = 0.9, 
n = 0.8, (e) A = 0.85, SI = 0.8, and (f) A = 0.8, O = 0.8. 
The darker color area indicates the lower density. 



which is similar to the case of SI = 0.59 ll27l . By com- 
parison, the l z in the equilibrium state for the former case 
is larger than that for the latter case (see Fig. 3 in ll27l "). 
which is caused by the further approach of the visible vor- 
tices to the center of the DW potential and the increase 
as well as the further concentration of the hidden vortices 
distributing along the central barrier (the phase profile for 
Q = 0.65 is not shown here). With the increase of SI, 
the visible vortices enter the BEC with the pair shape and 
form linear vortex lattice pair as shown in Figs. 4(b)-4(d). 
When SI = 0.85, the mirror symmetry is broken and a 
zigzag vortex lattice pair is formed [see Fig. 4(e)]. Finally, 
the transition into a triangular vortex lattice pair occurs, 
which can be seen clearly in Fig. 4(f) and Fig. 1(f). Here a 
common characteristic is that these visible vortex patterns 
exhibit the centrosymmetry due to the even parity of the 
DW potential with respect to the spatial coordinate. 

In Fig. 5, we present the steady density distributions 
\ip\ 2 at t — 250, where c = 600 and A takes on different 
values. For the fixed angular frequency SI = 0.9, when 
A = 2 the visible vortices form a linear lattice pair along 
the x axis due to the tight confinement in the y-direction of 
the trap as shown in Fig. 5(a). When A = 1.5, the visible 
vortices become a zigzag vortex lattice pair [see Fig. 5(b)]. 
A further reduction of A leads to a formation of a triangu- 
lar lattice pair of the visible vortices [see Fig. 5(c) and Fig. 
1(f)]. The visible vortex patterns for the fixed angular fre- 
quency SI — 0.8 with different values of A are illustrated 
in Fig. 4(d) and Figs. 5(d)-(f). With the decrease of A, the 
visible vortices entering the BEC by means of vortex pair 
arrange in a linear lattice pair along the y-direction be- 
cause of the strong restraint in the ^-direction of the trap. 
In particular, the two halves of the BEC approximate two 



rectangles at A = 0.8, where the linear visible vortex lat- 
tice pair are similar to two ionic crystals [45]. 

From Fig. 3(f) in Ref. l27l and Fig. 4(a), we can see 
that for 0.59 < SI < 0.65 there are always two visible vor- 
tices, i.e., the number of the visible vortices is degenerate 
with respect to the rotation frequency. Similar degeneracy 
also occurs in other certain intervals of rotation frequency. 
In fact, we find that the visible vortex number is possi- 
bly degenerate concerning to not only Q but also A, for 
instance, there are all the twelve visible vortices for the 
cases of 0.8 < A < 0.85 with the same SI = 0.8 as shown 
in Figs. 5(e)-5(f). The vortex patterns in Fig. 4 and Fig. 
5 are rather stable as verified in our simulation. Recently 
the stability of linear vortex lattice in a rotating anisotropic 
harmonic trap has been discussed in Ref. |[T6l . 

The structure change of vortex patterns in a BEC in a 
rotating anisotropic harmonic trap has been studied by sev- 
eral theoretical groups |Hl[l5][T7][l8][l9]g0]EI]E21, where 
the authors adopted mainly the LLL approximation or 
Thomas-Fermi (TF) approximation to analyze the vortex 
structures. However, all of these works did not consider 
the inevitable dissipation of rotating BECs. In addition, 
the analytical approaches in Refs. lfT^n3irT7irT8irT9l[20l[2Tl 
1221 can not be used directly to investigate the steady-state 
structures of a rotating DW BEC due to the presence of 
the Gaussian barrier. Here, for the rotating DW BEC there 
are three different kinds of vortices: visible vortices, hid- 
den vortices and ghost vortices, and each of them have 
their own unique topological structures [see Fig. 1(f), Fig. 
2(f), and Fig. 4], which is very different from the cases of 
rotating anisotropic harmonic trap. Moreover, as the val- 
ues of all the parameters are provided above, the structural 
change of the visible vortex patterns for the rotating DW 
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BEC with different Q or different A is experimentally ac- 
cessible and examinable. 

On the other hand, the linear visible vortex lattice pairs 
in Fig. 4 and Fig. 5 have potential applications in quan- 
tum information as they are somewhat similar to the ionic 
crystals l45l . For instance, a visible vortex in linear ar- 
rangement may be encoded as a quantum bit via a vortex- 
antivortex superposed state |j37]Sn]|46]|47]|48] . The flip be- 
tween the vortex and the antivortex can be realized either 
by optical methods 14911501 or by compressing and reha- 
bilitating the trapping potential l48l . The coherence time 
of vortex quantum bit is long enough to allow a series of 
operations lfl6l . In view of the degeneracy of visible vor- 
tex number, the small errors of fl or A due to experimental 
uncertainties do not influence the manipulation toward the 
vortex bits, which is beneficial to the information process- 
ing. Another advantage of the linear visible vortex lattices 
is that each visible vortex has only two nearest neighbors 
with controllable interaction, where the effect caused by 
the distant hidden or ghost vortices can be ignored because 
the particle density in the regions of the potential barrier 
and the outskirts of the cloud is very small. 

5. Conclusion 

In summary, we have numerically studied the dynamics 
of vortex formation and the structure of vortex lattices 
in a rotating DW BEC by using a phenomenological 
dissipation model 1271 . The vortex formation process in 
a rotating DW BEC is evidently different from that in a 
rotating single-well BEC 11261 . For the latter case, there 
are at most two types of vortices: visible vortex and ghost 
vortex. For our case, there exist three kinds of vortices: 
visible vortex, ghost vortex and hidden vortex. Once 
the DW potential rotates, the BEC undergoes complex 
turbulent oscillation principally consisting of the surface 
mode with I = 4 instead of the usual elliptic oscillation 
consisting of quadrupole mode with I = 2 in a rotating 
harmonic trap 11261 , The hidden vortices and ghost vortices 
always form earlier than the visible vortices. The critical 
frequency for the visible vortex generation is marked by 
a sudden jump in the evolution of the average angular 
momentum per atom and an inflexion in the evolution of 
the average energy per atom. The steady vortex patterns 
of the rotating DW BEC display well centrosymmetry 
because of the parity effect. It is shown that a structure 
change of the visible vortex patterns (including the 
formation of linear and zigzag vortex lattice pairs) can 
be realized by governing the rotation frequency and the 
anisotropic parameter of the DW potential. The structural 
transition of visible vortex patterns makes it possible 
to artificially manipulate the distribution of vorticity in 
superfluids and apply the quantized vortices to quantum 
information processing. 
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